//Only for re-run		
clear all
use "[Folder]\Panel_Firm_Lab_DirVarInd_temp.dta"

merge 1:1 LopNr_PeOrgNr year using "[Folder]\IndFirm_nr.dta", keepusing(Ind4_nr_q )
drop if _merge==2
drop _merge



	global MLab = 5
	global FLab = 2
	distinct Ind if u==0
	global Inds = r(ndistinct) - 1
	global Sizes = 4
	global Cities = 42
	global Years = 6
	global Res = 13
	global NbrOfParam =  $MLab + $FLab + $Cities + $Sizes + $Years + $Inds + $Res
	disp $NbrOfParam


//Function evaluation program
     program nlceslabrho
                version 14.2
                syntax varlist [if], at(name)
                local log_VA: word 1 of `varlist'
				
	// Retrieve parameters out of at matrix				
				tempvar LabAgg1 LabAgg2  temp_VA
				gen double `LabAgg1' = 0 `if'
				gen double `LabAgg2' = 0 `if'
                gen double `temp_VA' = 0 `if'
				
				replace `LabAgg1' = `LabAgg1' + RefGr  `if'

				forvalues i=1(1)2 {
				tempname delta 
				scalar `delta' = `at'[1,`i']
				replace `LabAgg1' = `LabAgg1' + `delta' * (lab_M_`i')  `if'
				}
				
				forvalues i=1(1)3 {
				tempname delta 
				scalar `delta' = `at'[1,2+`i']
				replace `LabAgg2' = `LabAgg2' + `delta' * (lab_M_C`i')   `if'
				}

				
				tempname  delta 
				scalar `delta' = `at'[1,$MLab + 1]
				replace `LabAgg1' = `LabAgg1' + `delta' * (lab_FNat) `if'

				tempname  delta 
				scalar `delta' = `at'[1,$MLab + 2]
				replace `LabAgg1' = `LabAgg1' + `delta' * (lab_FFor)  `if'
	
				forvalues i=1(1)$Inds {
				tempname Industry 
				scalar `Industry' = `at'[1,$MLab + $FLab + `i']
				replace `temp_VA' = `temp_VA' + `Industry'*Ind_`i' 	`if'			
				}

				//Adj index for number of firm size dummies
				forvalues i=1(1)$Sizes {
				tempname Size 
				scalar `Size' = `at'[1, $MLab + $FLab  + $Inds +`i']
				replace `temp_VA' = `temp_VA' + `Size'*FirmSize_`i' `if'				
				}
				//Adj index for number of City dummies
				forvalues i=1(1)$Cities {
				tempname City 
				scalar `City' = `at'[1, $MLab + $FLab  + $Inds + $Sizes +`i']
				replace `temp_VA' = `temp_VA' + `City'*k_`i' 	`if'			
				}
				//Adj index for number of year dummies
				forvalues i=1(1)$Years {
				tempname Year 
				scalar `Year' = `at'[1, $MLab + $FLab  + $Inds + $Sizes +$Cities +`i']
				replace `temp_VA' = `temp_VA' + `Year'*y_`i' `if'
				}				
				
	
                tempname  thetaL Cap Cap_Sq Cap_Qu TI TI_Sq TI_Qu TICap TICap_Sq TI_SqCap TI_SqCap_Sq ic  sigma

				scalar `thetaL' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +1]
				scalar `Cap' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +2]
				scalar `Cap_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +3]
				scalar `Cap_Qu' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +4]
				
				scalar `TI' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +5]
				scalar `TI_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +6]
				scalar `TI_Qu' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +7]
				
				scalar `TICap' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +8]
				scalar `TICap_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +9]
				scalar `TI_SqCap' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +10]
				scalar `TI_SqCap_Sq' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +11]

		
				scalar `ic' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +12]		
				scalar `sigma' = `at'[1,$MLab + $FLab  + $Inds + $Sizes +$Cities + $Years +13]		
				
				replace `temp_VA' = `temp_VA' + `ic'  + `thetaL' *(1/`sigma') *ln( `LabAgg1'^(`sigma') + `LabAgg2'^`sigma' )  + `Cap'*logCap + `Cap_Sq'*logCap_Sq + `Cap_Qu'*logCap_Qu + `TI'*logGrosInvTot + `TI_Sq'*logGrosInvTot_Sq  + `TI_Qu'*logGrosInvTot_Qu + `TICap'*logGrosInvTot*logCap + `TICap_Sq'*logGrosInvTot*logCap_Sq  + `TI_SqCap'*logGrosInvTot_Sq*logCap + `TI_SqCap_Sq'*logGrosInvTot_Sq*logCap_Sq  `if'
				replace `log_VA' = `temp_VA' `if'
            end
			///////////////////PROGRAM END///////////////////////////
			
			
 


			
//run regression, with number of param matching M			
mat M = J(1,$NbrOfParam ,1)

mat M[1,$NbrOfParam ]=0.7



//All firms
nl ceslabrho @ logFtg_Foradlingsvarde if u==0, nparameters($NbrOfParam ) initial(M) iter(75) variables(lab_M_C1-lab_M_C3 lab_M_1 lab_M_2 logCap logCap_Sq logCap_Qu logGrosInvTot logGrosInvTot_Sq  logGrosInvTot_Qu logGrosInvTot y_1-y_6 k_1-k_41 FirmSize_1-FirmSize_4  Ind_* lab_FFor lab_FNat  RefGr)
margins, dydx(RefGr lab_M_C*) atmeans
////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////
//Parameter values in table 8 are marginal effects for respective labor type normalized by the marginal effect of the reference group (except for sigma which is the parameter value from the original estimation)
////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////

//Upper NR-half
nl ceslabrho @ logFtg_Foradlingsvarde if u==0 & (Ind4_nr_q==3 | Ind4_nr_q==4) , nparameters($NbrOfParam ) initial(M)  iter(75) variables(lab_M_C1-lab_M_C3 lab_M_1 lab_M_2 logCap logCap_Sq logCap_Qu logGrosInvTot logGrosInvTot_Sq  logGrosInvTot_Qu logGrosInvTot y_1-y_6 k_1-k_41 FirmSize_1-FirmSize_4   Ind_* lab_FFor lab_FNat  RefGr)
margins, dydx(RefGr lab_M_C*) atmeans
////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////
//Parameter values in table 8 are marginal effects for respective labor type normalized by the marginal effect of the reference group (except for sigma which is the parameter value from the original estimation)
////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////

//Lower NR-half
nl ceslabrho @ logFtg_Foradlingsvarde if u==0 &(Ind4_nr_q==1 | Ind4_nr_q==2), nparameters($NbrOfParam ) initial(M) iter(75) variables(lab_M_C1-lab_M_C3 lab_M_1 lab_M_2 logCap logCap_Sq logCap_Qu logGrosInvTot logGrosInvTot_Sq  logGrosInvTot_Qu logGrosInvTot y_1-y_6 k_1-k_41 FirmSize_1-FirmSize_4   Ind_* lab_FFor lab_FNat  RefGr)

margins, dydx(RefGr lab_M_C*) atmeans
////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////
//Parameter values in table 8 are marginal effects for respective labor type normalized by the marginal effect of the reference group (except for sigma which is the parameter value from the original estimation)
////////////////////////////////////////////////////////////////////////////////////
////////////////////////////////////////////////////////////////////////////////////
